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<^ ; Abstract 

■ Recent progress in unitarity techniques for one-loop scattering amplitudes makes a numerical 

1—1. implementation of this method possible. We present a 4-dimensional unitarity method for calculat- 

Qh| ing the cut-constructible part of amplitudes and implement the method in a numerical procedure. 

g_) ' Our technique can be applied to any one-loop scattering amplitude and offers the possibility that 

^ : 

one-loop calculations can be performed in an automatic fashion, as tree-level amplitudes are cur- 

, rently done. Instead of individual Feynman diagrams, the ingredients for our one- loop evaluation 

OO ' 

Q>^ \ are tree-level amplitudes, which are often already known. To study the practicality of this method 

m ■ 

^SJ ■ we evaluate the cut-constructible part of the 4, 5 and 6 gluon one-loop amplitudes numerically. 



using the analytically known 4, 5 and 6 gluon tree-level amplitudes. Comparisons with analytic 



OO 

O 

, 

, answers are performed to ascertain the numerical accuracy of the method. 
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I. INTRODUCTION 



Analytic unitarity techniques in Feynman diagram calculations have been used for a long 
time [l|,l2,y]- Their use in the context of gauge theories is even more powerful 0, Q] and they 
were successfully applied to the calculation of one-loop amplitudes of phenomenologically 
important 5-leg and 6-leg processes in QCD (for a recent review see I?]). 

In gauge theories the conventional Feynman diagram method produces intermediate re- 
sults which are much more complicated then the final answer. One evaluates the numerous 
non gauge-invariant individual Feynman diagrams by expanding the tensor loop integrals 
into form factors. This decomposition generates a large number of terms. With a grow- 
ing number of external particles it becomes a forbidding task to simplify the expression 
analytically. This forces one to adopt more numerical techniques (see e.g. ref. js]), which 
can be computationally intensive due to the large number of terms. In addition, the large 
cancellations between the Feynman diagrams can potentially lead to numerical instabilities. 

As an alternative, the unitarity cut method uses only on-shell states, manipulates gauge 
invariant amplitudes and has been used to derive simple answers with simple intermediate 
steps {g]. New ideas on twistors 9], multipole cuts (generalized unitarity) 10 1, recursion 



relations 
dimension 
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17 



14|. algebraic reduction of tensor integrals [15|] and unitarity in D- 
191] have made the unitarity cut method even more promising. It 



appears that ultimately one can find an efficient algorithm which can be used to calculate 
the one-loop amplitudes in terms of tree-level amplitudes. The progress is due to three 
important observations. 

First, any one-loop amplitude can be decomposed in terms of scalar box, triangle and 
bubbl 
shell [20 

may have infrared (boxes and triangles) or ultraviolet divergences (bubbles). 

The second key observation concerns the application of unitarity techniques to amplitudes 
with multiple cuts. Using unitarity techniques the coefficients of the master integrals are 
determined by multiple cuts of the amplitude which place the cut internal lines on their mass 
shell. After cutting, the tree-level 3-gluon scattering amplitude with all 3 gluons on-shell 



e master integrals where both internal and external particles can be massive or off- 
2l[ |. The master integrals^ have to be calculated in dimensional regularization and 



^ For a collection of currently known one- loop master integrals, see the web-site http : / / qcdloop . f nal . gov 
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can appear at a vertex of the diagram. These 3-particle amphtudes are identically zero by 
momentum conservation. Therefore coefficients of the associated master integral can not be 
extracted. This obstacle is removed by the observation that the tree-level helicity amplitudes 
can be analytically continued to complex momentum values [o], 1^ allowing for solutions to 
the unitarity constraints in terms of the complex loop momentum. All of the relevant tree- 
level expressions are non-zero and the appropriate one-loop amplitude is reconstructible from 
the tree-level amplitudes. With this method the coefficients of the 4-point master integrals 
can be extracted both analytically and numerically. They are given in terms of the product 
of 4 tree- level amplitudes, evaluated with the complex on-shell loop momenta [l^ . However 
the coefficients of the 3- and 2-point master integrals were still difficult to extract because 
terms already included in the 4-point contributions had to be subtracted. It was not clear 
how to express this subtraction in terms of the corresponding tree-level amplitudes. 

The third important observation is that there is a systematic way of calculating the 
subtraction terms at the integrand level. By manipulating the one-loop amplitude before 
the loop integration is carried out, the unitarity method is reduced to the algebraic problem 
of a multi-pole expansion of a rational function. Alternatively, the method of ref. [l^ can 
also be viewed as the calculation of the residues of each pole term of the integrand. The 
resulting 4-propagator pole (i.e. box contribution), 3-propagator pole (triangle contribution) 
and 2-propagator pole (bubble contribution) naturally decompose the loop momentum inte- 
gration vector into a "physical" space spanned by the respective external momenta and the 
remaining "trivial" space orthogonal to the "physical" space. The so-called spurious terms 
(or subtraction terms) of ref. 15| are determined by the most generally allowed dependence 
of the residue on the components of the loop momentum in the "trivial" space. By definition, 
these spurious terms vanish upon integration over the loop momentum. These ideas allow 
the extraction of all the coefficients of the master integrals for a given one-loop amplitude. 
A possible algorithm for analytical extraction of the coefficients within the unitarity method 
has been worked out in ref. 22 1. 



In this paper we will expand on the algebraic method [15|] by developing a numeri- 
cal scheme. With the numerical method outlined in this paper we evaluate only the cut- 
constructible part of the amplitude. We will show that the master integral coefficients are 
calculated in terms of tree-level amplitudes. This makes it possible to "upgrade" existing 
leading order generators to produce the cut-constructible part of the one-loop amplitudes. 
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The "upgrade" requires allowing two of the external momenta in the tree-level amplitude 
to be complex 4-vectors, while leaving the analytic expression of the tree-level amplitude 
unchanged. For example, to evaluate the cut-constructible part of the 6-gluon amplitude 
only the analytic 3-, 4-, 5- and 6-gluon tree-level amplitudes are needed (expressed in spinor 
product language or any other form). Alternatively, one could use an efficient recursive 
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numerical method to evaluate the tree- level amplitude fll . 

Because we implement 4-dimensional unitarity cuts, the so-called rational part is not 
generated . In principle we could expand our scheme to D-dimensional unitarity cuts, 
thereby generating the complete amplitude. However, it is not immediately clear how to do 
this while maintaining the requirement of 4-dimensional tree-level building blocks. Alterna- 
tive methods exist to determine the rational part, which in principle can be combined with 
the numerical method outlined in this paper to give the complete scattering amplitude. A 
method for determining the rational part using on-shell recursion relation has been success- 
fully developed and used 2J, |25[. A direct numerical implementation of this method should 
be possible and especially attractive due to the recursive nature of this method. Direct 

nnnn 

numerical implementation of the D-dimensional unitarity method |16l . 1171 . 1181 . 119I | is harder 
but may also become practical in the near future. Other methods have been developed using 
Feynman diagram expansions. While these methods are attractive as far as simplicity goes, 



they re-introduce all the problems of standard Feynman diagram calculations. In ref. [26 1 
a method for calculating the rational part is developed partially based on Feynman dia- 
gram calculations. In refs. {27, 28 1 a method is proposed to use simplifed Feynman diagram 
techniques for calculating the rational part. 

In section 2 we will derive the formalism, which we then apply in section 3 to calculate 
the cut-constructible part of the 4-, 5- and 6-gluon amplitudes at one-loop by only using 
the tree-level amplitudes for the 4-, 5- and 6-gluon amplitudes. We will then compare the 
numerical results to the analytic calculations. 
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FIG. 1: The generic A^-point loop amplitude. 
II. THE STRUCTURE OF THE ONE-LOOP INTEGRAND FUNCTION 

The generic D-dimensional A^-particle one-loop amplitude (fig. [1]) is given by^ 

An{puP2,---,Pn) = J [dl\ ^^^^ , (1) 

where pt represent the momenta flowing into the amplitude, and [dl] = d^l. The numerator 
structure J^{pi,P2, ■ ■ ■ ,Pn', is generated by the particle content and is a function of the 
inflow momenta and the loop momentum. Since the whole amplitude has been put on a 
common denominator, the numerator can also include some propagator factors. The depen- 
dence of the amplitude on other quantum numbers has been suppressed. The denominator 
is a product of inverse propagators 

di = di{l) = {l + q,f -m^i = (^-qo + Y,p}j - , (2) 

where the 4- vector go represents the arbitrary parameterization choice of loop momentum. 
The one-loop amplitude in D = 4 — 2e can be decomposed in the scalar master integral basis 

^ Wc restrict our discussion to (color) ordered external legs. The extension for more general cases is 
straightforward . 
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giving ^ 



{PUP2, ...,Pn) 

l<jl<i2<«3<«4<A^ 



l^2^3 

l<ii<i2<i3<N 



+ X Ki2iPl^P2,---,PN)Ihi2 

l<ii<i2<N 

+ X «ii(Pl,1^2,---,P7v)/n , (3) 



l<n<Ar 

where the master integrals are given by 

4...M= [[dl] , \ . (4) 

Analytic expressions for the master integrals with massless internal lines are reported in 
ref. |5|. 

The maximum number of master integrals is determined by the dimensionality, D, of 
space-time; for the physical case this gives up to 4-point master integrals. The unitarity 
cut method is based on the study of the analytic structure of the one-loop amplitude. The 
coefficients are rational functions of the kinematical variables and will in general depend 
on the dimensional regulator variable e = (4 — D)/2. When all the coefficients of the 
master integrals are calculated in 4 dimensions we obtain the "cut-constructible" part of 
the amplitude. The remaining "rational part" is generated by the omitted 0{e) part of the 
master integral coefficients [7|. 

For a numerical procedure we need to recast the study of the analytic properties of the 
unitarity cut amplitudes into an algebraic algorithm which can be implemented numerically. 



In ref. 



15| it was proposed that one focus on the integrand of the one-loop amplitude, 



An{Pi,P2,---,Pn\1) = j-i 1 • (5) 

aia2 ■ ■ ■ cln 

This is a rational function of the loop momentum. Any A^-point tensor integral of rank M 
(M < A^) with A^ > 5 can be reduced to 4-point tensor integrals of rank K {K < 4) by 
application of Schouten identities. Therefore we can re-express the rational function in an 



^ We drop the finite 6-dimensional 5-point master integral because its coefficient is of ©(e) |29| and therefore 
it will not contribute to the final answer where we take e — > 0. 
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expansion over 4-, 3-, 2- and 1-propagator pole terms. The residues of these pole terms 
contain the master integral coefficients as well as structures which reside in the subspace 
orthogonal to the subspace spanned by the external momenta. These spurious terms are 
important as subtraction terms in the determination of lower multiplicity poles. The number 
of spurious structures is 1 for the box, 8 for the triangle, 6 for the bubble, and 4 for the 
tadpole. After integration over the loop momenta, Eq. ([3]) is recovered. This approach 
transforms the analytic unitarity method into the algebraic problem of partial fractioning 
a multi-pole rational function. The remaining integrals after the partial fractioning are 
guaranteed to be the master integrals of Eq. (jll). This makes a numerical implementation 
feasible. 



A. The van Neerven-Vermaseren basis 

Consider a set of R inflow momenta, /ci, . . . , /cj? in a D-dimensional space-time^. Taking 
momentum conservation into account, Xl^i = 0? the physical space spanned by the 
momenta ki has dimension min(Z},i? — 1). As a consequence, for R > D + 1 additional 
Schouten identities exist, which can be exploited to prove the D = 4 master integral basis 
of Eq. ([3]) 2l|. For R < D, the physical space forms a lower dimensional subspace. We can 



define an orthonormal basis, the van Neerven-Vermaseren (NV) basis 30|, which separates 
the D-dimensional space into the Dp-dimensional "physical" space and the orthogonal Dj- 
dimensional "trivial" space where 

D = Dp + Dr, Dp = min {D, R-1); Dt = max (0, - i? + 1) . (6) 



The inflow momenta are either equal to the external momenta pi, or to sums of external momenta. 
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To define the NV-basis we introduce the generalized Kronecker delta 31 1 



i:^^l^J.2■■■^J.R 



yi 1^2 



the compact notation 



"vi "i/2 • • • "ur 



^p^^2■■■^lR — s^ll^l2■■■^lRl. „V2 

'^uiq-'-UR — '^yxU2 - UR '^'MiV 



(7) 



and the {R — l)-particle Gram determinant 



(8) 



(9) 



Note that for R > D + 1 the generalized Kronecker delta is zero. For the special case 
D = R we have the factorization of the Kronecker delta into a product of Levi-Civita 
lensoib. ~ t ^iyiU2 -UR- 

Some examples of the generalized Kronecker delta are 



xkik2 



919293 



h ■ qi k2f, - k2-qi ki^ 

^i-^i^9'2t-^i-^2e + A:i-g3e 
+ki ■ qi {k2 ■q2k3-q3-k2- qs k^ ■ ga) 

-ki ■ q2 {k2 ■qik3-q3-k2- qs k^ ■ qi) 

+ki ■ qs {k2 ■ qik^- q2- k2- q2 ks ■ qi) ■ 



(10) 



We now want to construct the NV-basis for R momenta. We define Dp basis vectors 

^ki...ki-iiiki+i...kr,p 
u, / 1 1 \ ki...ki—±kiki^i...kjjp 



(11) 



^ This notation is closely related to the asymmetric Gram determinant notation of ref. 32 1 

fci • • • kji 



G 



51 • • • QR 
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with the properties Vi ■ kj = Sij for j < Dp. When R < D we also need to define the 
projection operator onto the trivial space 

A{ki, . . . , kR-i) 

with the properties w^'^ = Dt = D + 1 — R, k^ w^u = and w^aW°'^ = w^" . Note that this 
operator is the metric tensor of the trivial subspace, with the decomposition 

w^^= J2 (13) 

i=l 

where the D + 1 — R orthonormal base vectors of the trivial space rij have the property 
rii ■ Uj = 6ij, rii ■ kj = rti ■ Vj = 0. 

The full metric tensor decomposition in the NV-basis is given by ^ 

Dp Dp Dt 

= E + ^"'^ = E + E « • (14) 

i=l i=l i=l 

For the case D = R the sole basis vector of the 1-dimensional trivial space is proportional 
to the Levi-Civita tensor. For the cases R < D we can explicitly construct the basis vectors 
fulfilling all the requirements. 

As an example in the case of D = 4 and i? = 4 we get 

srfik2k3 Sjkink^ diikiH 

k,, k,) = <(^i' ^2, k,) = ^3 (^1, h, k,) - '^'^'^ 



Aik,,k2,ksy ^' A{k,MMy " ' HkiM.h) 

^kik2k3U kik2k3U 

We want to decompose the loop momentum into the NV-basis for a graph with denom- 
inator factor di, (^2, . . . , d^. The denominators are as usual given by di = (/ + g,/)^ — 
and ki = Qi — By contracting in the loop momentum with the metric tensor given in 

Eq. f[T^ we get the loop momentum decomposition in the NV-basis 

Dp Dt 

F = J2l-kiv'^ + Y,l-n,n'^ . (16) 

1=1 i=l 

Using the notation I ■ Ui = ««(/) = ai and the identity 

l-ki = ^[di- di^i - (g2 - m^) + _ ^2_jj ^ ^^^^^ 



^ By expanding the generalized Kronecker delta functions in the Vi veetors one can show that k'j^Vi 
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we find 

I' = ^« + E ^(^^ - ^-i) < + E «^ < ' (18) 

i=l i=l 

where do = dn, rriQ = rriR and 

i=l 

In the case that R > D + 1 the decomposition of the loop momentum into the NV-basis 
implicitly proves Eq. ([3]). Also, when R < D it allows us to include the unitarity constraints 
without resorting to the explicit 4-dimensional spinor formalism used in analytic calculations. 
By avoiding the 4-dimensional spinor formalism, the formulation is also valid for massive 
internal particles (where the mass can be real or complex valued). 

For example, in the case of a 4-dimensional pentagon, [D = 4 and i? = 5) we get 

P = Vi" + l{di- d,) < + - d^) v!^ + ^(4 - ^^2) + ^{d,- ds) < 

T//i 1/ 2 2 2 I 2\ U ^ /■ 2 2 2 I 2\ U 

H = -gi^i + "^s) < - 2(^2 - gi - ^2 + mj vi; 

- ?2 - "^3 + T^l) V3 - ^{ql -ql-ml + ml) < . (20) 

Similarly for a 4-dimensional triangle {D = A and i? = 3) we get 

P = V^3^ + ^(rfi-rf3)< + ^(rf2-c?i)^^2^ + ai< + a2ri^ 
^3" = -^('?? - ql + ml) < - \{ql - - + m?) ^;2^ . (21) 

Thus we see that the same basis decomposition is used for the tensor reductions in the case 
R> D + 1 and solving the unitarity constraint in the case that R < D. 



B. Partial fractioning of the integrand 

For the remainder of the paper we restrict ourselves to a 4-dimensional space. Given 
the master integral decomposition of Eq. ([3]) we can partial fraction the integrand of any 
4-dimensional A^-particle amplitude as 

i<h<i2<h<u<N diA2di3du i<ij<i2<i3<jv diidiAs i<i^<i^<N ^'^^'^ i<h<N 
To calculate the numerator factors, we will calculate the residues by taking the inverse 
propagators equal to zero. The residue has to be taken by constructing the loop momentum 
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lij...k such that di{lij...k) = dj{lij...k) = ■■■ = dk{lij...k) = 0. Then the residue of a function 
F[l) is given by 

ReSij...k[F{l)]= U{l)d,{l)---dk{l)F{l)) _ . (22) 
The specific residues are now given by 
dijki{l) = Resijfczf^7v(0) 

dijki{l) 



Cijk{l) = Resijk An{1) - ^ 



, , . . , didjdkdi 

lfi,j,k ■' 

Cijk{l) 1 \ ^ dijki{l) 



fe.,(/) = Res.,U.(0-ETf^-^E 



didjdk 2! /-f, . didjdkdi 

K^i,] ■' k,lfi,] 

bij{l) 1 Cijk{l) 1 \- dijki{l) 



,0. Re, (4.(0-5:^-^5:^-^ E 



didj 2! ^.didjdk 3! .f-^. didjdkdi 



(23) 



Note that the coefficients are defined to be symmetric in the propagator indices (e.g. C123 = 
C213 = C312) and coefficients with repeated indices are to be set to zero (e.g. (iiase = C112 = 0). 
As an example, some residues of a 5-particle amphtude are given by 



dl245 = ReSi245 (^^Af(0) 



-D I A n\ ^1235(0 C?2345(0 

C235 = ReS235 An{1) - 



bi4 = Resi4 ( An{1) 



did2d^d5 d2d^did5 

Cl24(0 Ci34(/) Ci45(/) ^1234(0 ^^1245(0 '^1345(0 



did2d4 did^di did^d^ did2d-idi did2di^d^ did^d^d^ 



(24) 



In the following sub-sections we will explicitly construct the residue functions using only 
tree-level amplitudes. This construction is well-suited for numerical implementation. 



C. Constructing the box residue 

To calculate the box coefficients we choose the loop momentum lijki such that four inverse 
propagators are equal to zero, 

dijkiihjki) = RGSijki(^AN{l)^ ■ (25) 
11 



We will drop the subscripts on the loop momentum in the following. Because we have to 
solve the unitarity constraints explicitly, we have to choose a specific parameterization, qo, 
in Eq. ([2]). Using the NV-basis of the four infiow momenta for the box and using the fact 
that di = dj = dk = di = we can use Eq. (ITSl) to decompose the loop momentum as 

//^ = + ai < . (26) 

Choosing for the parameterization go = Yl]=iPj (^^ usual the index i is understood to be 
modulo A^) such that = we have 

^4 = -^1 + < - ^il'j - li - + "^?) ^2 ^1 - "^fc + "^?) ^3 > (27) 
where 

r/ifeafca rfclM^S rfclfc2M ^fj.kik2k3 

M. _ fclfc2fc3 At _ fcifc2fc3 ^^ _ kik2k3 a _ ^ /r)o\ 

and 

ki = Qi-, k2 = qj - qi] ks = qk- qf, A{K k2, k^) = S^l'^f^l . (29) 

The variable cti will be determined such that the unitarity condition d^ = dj = dk = di = 
is fulfilled^. Imposing the constraint dn = for n = {i,j, k, 1} using Eq. fl26l) 

dn = = (/ + g„)^ - ml 

= (1^4 + aiUi + qnf - ml 

= Vl + al + 2V^ -qn + qn-ml, (30) 



and using the identity 



so that 



Vi ■ qn = -^{ql - ml + m^) , (31) 



al = -{V,' - mf) , (32) 



we find two complex solutions 



l^ = Vr±iJVi-mfxnt, (33) 



^ In fact ai = cti{l), by varying the loop momentum (allowing for complex values) wc can set ai to the 
required complex value. This dependence is implicitly assumed, we will treat ai as a "free variable" . 
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which are easily numerically implemented. We note that because di = dj = dk = di = 0, we 
have (/J)^ = ml for n = z, j, k, I where = if-j^^ + Qn = ifjki + YJj=i+iPi- other words 
the four propagators are on-shell and the amplitude factorizes for a given intermediate state 
into 4 tree- level amplitudes AA'^'^\ For the residue of the amplitude in Eq. fl2^ we find (all 
indices are assumed modulo A^, i.e. i = n + N = n) 

ReSi,ki[AN{l^)) = M^'\lt,P^+l, . . .,Pf, -If) X M^'\lf;p,+u . . . 

X M^'\lt,Pk+u ■■■,Pi; -If) X M'^'\lf;pi+i, . . . -if) , (34) 

where the loop momenta are complex on-shell momenta and there is an implicit sum 
over all states of the cut lines (such as e.g. particle type, color, helicity). For example, the 
residue of the amplitude for the pure 6-gluon ordered amplitude with dg = d2 = d^ = d^ = 
factorizes into (see fig. [2]) 

Res2346 (A(/^)) = Mf\li;puP2;-lt) y< Mf\lt,ps;-li) 

X Mf\lf;p,; -It) X MT{lhp,,P,\ -if) , (35) 

where the implicit sum over the two helicity states of the four cut gluons is assumed. The 
tree-level 3-gluon amplitudes, M.f \ are non-zero because the two cut gluons have complex 
momenta |10[ |. 

Any remaining dependence of the residue dijki on the loop momentum enters through its 
component in the trivial space, 

'dijki{l) = dijki{ni ■ I) . (36) 

The number of powers of the loop momentum / in the numerator structure is called the rank 
of the integral. After integration we find using Eq. ( fT3l) that {rii ■ /)^ ~ = 1. Thus rank 
one is the maximum rank of a spurious term (which by definition vanishes upon integration 
over /). Hence the most general form of the residue is 

dijki{l) = dijki + dijki I ■ rii . (37) 



Using the two solutions of the unitarity constraint, Eq. (l33l) . we now can determine the two 
coefficients of the residue 

ReSijki(^AN{l^)^ +ReSijki(^AN{l~ 



dijki 



2 



Resijki (^An{1^)^ - Resijki (^An{1'^ . ^ 
^ijki = ^. /-.^o =s ■ (38) 
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4 



P3^ 



P4 



l3=l+Pl+P2+P3 



l2=l+Pl+P2 i 


1 s» 

I j 






1 

^ ^ ' 


' l4=l+Pl+P2+P3+P4 








Ps^ P/ 


i 


^ Pe ^ Ps 



FIG. 2: The factorization of the 6-gluon amphtude for the calculation of the ^2346 (0 residue with 
the loop momentum parametrization choice qq = 0. 

With the above prescription it is now easy to determine the spurious term for any value 
of the loop momentum. Finally we note that the integration over the term 



[dl] 



[dl] 



dijki + dijki Til ■ I 



dijki / [d I] 



dijkllijkl 1 



didjdj^di J didjdf^di J didjdf^di 

is now trivially done, giving us the coefficient of the box times the box master integral. 



(39) 



D. Construction of the triangle residue 

To calculate the triangle coefficients we need to put three propagators on-shell. Care 
has to be taken to remove the box contributions by explicit subtraction. Thus, the triangle 
coefficient is given by 



c,,,(/) = Res,,, I A^{1) - J2 fif ) • (40) 

l^i,j,k * 



Decomposing the loop momentum in the NV-basis of the three inflow momenta of the triangle 
with di = dj = dk = (choosing = 0) gives us according to Eq. (fTSl) 

//^ = + ai< + a2n^ , (41) 
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FIG. 3: The factorization of the 6-gluon ordered amphtude for the calculation of the 0234(1) residue 
with the loop momentum parametrization choice qq = —p^ — Pq = pi + P2 + Ps + Pa- 



with 



where 



and 



Vs' = -^iq^ -m^ + ml) < - - - m] + m^) vl^ , (42) 



r^fe skill 
%fc2 . Li ^ "fcifc2 

A{h,h)' ' A{h,k2) 



(43) 



h = qu h = qj - qu A{h, ^2) = • (^4) 

The base vectors of the trivial space {ni,n2} have to be explicitly constructed using the 
constraints 

Hi ■ Uj = 5ij] Hi ■ kj = 0; w^^{ki, k2, k^) = n^n^ + 77-2^2 • (45) 
The unitarity constraints [di = dj = dk = 0) give an infinite set of solutions * 

I'aia, = Vr + ai < + a2 n^; af + = -{V^ - ml) . (46) 



For masslcss internal lines the parameterization of ref. 2j| is obtained by taking ai = k x (a + i 6) and 
a2 ~ K X {a — ib) where = ^^3^1 o, = i(t — 1/t) and b = i(< + 1/t). By taking the t ^ 00 limit one 
gets the coefficient of the triangle master integrals. 
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In this case the residue of the amphtude factorizes into three tree-level amplitudes 

Res,,fc(^Ar(r"^)) = A^(°)(/r°^;Pm, • • ■,Pr, -If') x M^'\lf'''-p,^u • • • -1^') 
X , (47) 

with an implicit sum over the internal states of the cut lines. For example, the residue of 
the amplitude for the pure 6-gluon amplitude with = = ^4 = factorizes into (see 

fig. ED 

Res234(A(r^"^)) = -Mf (^2^"^P3;-^3'"^) xA^f (/3"^'^^P4;-C"') 

X Mf{lT';p,,p,,p,,P2-,-q^^') . (48) 

The remaining dependence of the residue Cijk on the loop momentum resides in the trivial 
space 

Cijk{l) = Cijki'Si, sa); si = rii ■ I, S2 = n2 ■ I . (49) 

The maximum rank of the triangle diagrams in standard model processes is 3. This gives 
us 10 possible terms ({1, si, S2, sf, siS2, si, sf, sfs2, sis^, s^}) for the most general polynomial 
form of the residue. However using Eq. fll3p we have the constraint + s| ~ + n| = 2 
which reduces the number of terms to 7. The form we chose is given by 

C^Al) = + Cg.i + CQS2 + C^lisl - sj) + .^.^(cg + C^S, + <^S2) . (50) 

We need to determine all 7 constants c-"^ by constructing 7 equations. This is accomplished 
by choosing 7 combinations of (ai, 0:2) with af + a^ = —{V-^ Eq. (1401) . This system 

of equations can be easily solved using a matrix inversion. Note that in principle we can 
generate an unlimited set of equations. This can be useful in a numerical application to 
obtain better numerical accuracy. 

With the above prescription it is now easy to determine the spurious term for any value 
of the loop momentum. Finally we note that the integration over the term 

is now trivially done, giving us the triangle coefficient times the triangle master integral. 
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E. Construction of the bubble residue 



This sub-section follows closely the previous two sub-sections. To calculate the bubble 
coefficients we need to put two propagators on-shell. The box and triangle contributions 
need to be explicitly subtracted. This gives for the bubble coefficient 



MO^Re=« ^^(O- Eg^-^ Eii^ . (62) 

Decomposing the loop momentum in the NV-basis of the two inflow momenta with di = 
dj = (choosing qj = 0) gives us according to Eq. fITS]) 

F = V2 + ai?7,^ + 02^2 + a3?T-3 , (53) 

with 

V2' = -l{ql-m^ + m'^)v^, , (54) 

where 

~ A{h) ~ kj ' ^^^^ 

and 

ki = Qi ■ (56) 

The base vectors of the trivial space {rii, n2, n^} have to be explicitly constructed using the 
constraints 

rii ■ rij = 5ij] rii ■ kj = 0; w^'^{ki, k2) = riiTil + 722^2 + n'^n'^ ■ (57) 

The solution to the unitarity constraints {di = dj = 0) gives as a infinite set of solutions 

^aiaaas = + «i < + «2 ni^ + u'^] aj + + = -{Vi - rr?j) . (58) 

As before the residue of the amplitude factorizes into tree-level amplitudes, in this case into 
2 tree-level amplitudes with an implicit sum over the states of the cut lines 

Res,,(^^(r^"^"^)) = 

-M(°nr"'"';Pm, • • • -/f "^"^) X A^W(/J^°^°^p,+i, . . . -/^i"2°3) _ (59) 

For example, the residue of the amplitude for the pure 6-gluon amplitude with d2 = d^ = ^ 
factorizes into (see fig. Hj) 

Res24(A(r^"^"^)) =A^fH/r"'"^P3,P4;-/r"^"^) X Alf (C"^"^P5,P6,Pl,P2;-/2^ . 

(60) 
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FIG. 4: The factorization of the 6-gluon amphtude for the calculation of the 624 (/) residue with 
the loop momentum parametrization choice qq = — ps — Pe = Pi + P2 + Ps + P4- 

The remaining loop dependence of the residue bij is in the trivial space 

bijil) = bijisi, S2, S3); si = ni ■ /, S2 = na ■ /, S3 = ng ■ / . (61) 

The maximum rank of the bubble in standard model processes is 2. This gives us 10 possible 
terms ({1, Si, S2, S3, s^, s^, S3, S1S2, S1S3, S2S3}) for the most general polynomial form of the 
residue. However, using Eq. (fT3|) we have the constraint sf + + si nl + + nl = 3 
which reduces the number of terms to 9. The form we chose is given by 

%(0 = C + ^Sj^^l + ^Sf ^2 + S3 + 6lf (S? - S^) + 6Sf (S^ - S^) + foSf S1S2 ^ 

(62) 

We need to determine all 9 constants fe^-"^ by constructing 9 equations. This is accomplished 
by choosing 9 combinations of (ai,Q;2,a3) with af + a^ + al = —(y2—m'^) in Eq. (l52l) . This 
system of equations can be easily solved using a matrix inversion. Note that in principle we 
can generate an unlimited set of equations. This can be useful in a numerical application 
to obtain better numerical accuracy. We only need the full bubble residue in the case when 
the tadpole contribution is non-zero. 

With the above prescription it is now easy to determine the spurious term for any value 
of the loop momentum. Finally we note that the integration over the term 

/<" If -'S'/^' 4^ 
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is now trivially done, giving us the bubble coefficient times the bubble master integral. 
F. Construction of the tadpole coefficient 

In the case there is a tadpole contribution we need to determine its coefficient from the 
relation 



a,; 



Using Eq. ( fTSj) the loop momentum can be decomposed in four orthonormal vectors 

= ttin^ + a2^2 + '^s'^s + tt4'^4 5 (65) 

and 

rii ■ rij = Sij] g'^" = n^n\ + n^n^ + n'^nX, + ^4^4 . (66) 
The solution to the unitarity constraints {di = with qi = 0) gives an infinite set of solutions 

2 I 2 I 2 I 2 2 ^aiy\ 

+ a2 + + = nij^ , (67) 

The residue of the amplitude becomes a tree-level amplitude with an implicit sum over the 
states of the cut particle 

ReSi(AN{l"'"'"'"')) = A^(°H^»°'°'°'°';Pi+i, • • • ,Pi; -/fi"2"^"*) . (68) 

For example, the residue of the amplitude for the pure 6-gluon amplitude with ^5 factorizes 
into (see fig. [5]) 

ReS5(A(r^"^"'^"^)) = -Mf (r"'"^"^P6,Pl,P2,P3,P4,P5; -/-^"^"3«4) _ (gg) 

The maximum rank of the tadpole in the Standard Model is one, giving for the spurious 
term 

a, = af ) + + a^s, + af^s, + a'fh, . (70) 
The coefficient of the master integral Jj is easily obtained by e.g. 

(71) 



(0) ai(ri°2a3a4) _^^.(^//3ia2a3"4) 



2 

with «! = — /?! = nii and ^2 = Q^s = ^4 = 0. 
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Pa' 



FIG. 5: The factorization of the 6-gluon amphtude for the calculation of the a^{l) residue with the 
loop momentum parametrization choice qo = —pQ = pi + P2 + P3 + Pa + P5- 

III. NUMERICAL RESULTS 



As an application we calculate the 4, 5 and 6 gluon scattering amplitudes at one-loop 
using the method of sec. 2. The cut-constructible parts of the ordered amplitudes are also 



known analytically ( 3J, l35|, l36|] , 37|] 



m 



38 



39 



40 



41 



42|). making a direct comparison 



possible. These multi-gluon scattering amplitudes form a good test for numerical procedures 
as they are the sum over a large number of Feynman graphs with significant gauge cancella- 
tions. Also, the 6-gluon amplitude was numerically evaluated using the integration-by-parts 
...odfr. Then— eva,„a*„W using that method was around 9 seconds per 
ordered amplitude (on a 2.8GHz Pentium processor). We can compare this evaluation time 
with the unitarity method of the previous section. This directly compares the computational 
effort between a numerical method using Feynman diagrams with form factor expansion and 
the numerical unitarity method using analytical expressions for the tree-level amplitudes. 



This Feynman diagram calculation yielded the full amplitude including rational terms. 
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To calculate the gluon scattering amplitude we need to determine all the master inte- 
gral coefficients and combine these with the master integrals according to Eq. ([3]). It is 
straightforward to numerically evaluate a coefficient using the method outlined in the pre- 
vious section. Given a set of external momenta we simply calculate the appropriate four 
vectors {vi} and {rii} according to Eqs. (1251 H5l H5| IFTI) . From these vectors we construct 
the special loop momenta of Eqs. (!26l 1411 [53l) which sets the appropriate denominators to 
zero. Using the analytically known leading order gluon amplitudes we can calculate the 
coefficient of Eqs. (138| [STl [63l) . Note that the multi-gluon scattering amplitudes do not get 
a contribution from the tadpole master integrals. We therefore need only to calculate the 
6^°^-coefficient of the bubble residue and not the remaining 8 coefficients of the spurious 
term. 

The ffist check on the numerical implementation is performed by calculating the e^^ term, 
where D = A — 2e. This term gets contributions both from the box master integrals and 
from triangle master integrals with one leg off-shell. The e~^-term for a n-gluon ordered 
amplitude is proportional to the leading-order amplitude and is given by 



m(i)(l,2,...,n)~-^xm(°)(l,2,...,n) + 0(e-^) . (72) 



The next check is the term. Again, the term is proportional to the leading-order 
ordered amplitude. This term gets contributions from the one-, two- and three-leg off-shell 
box master integrals, one- and two-leg off-shell triangle master integrals and bubble master 
integrals. The contribution is given by 



2, . . . , n) ~ (^-^ + ^ (^-y + log j j xmW (1,2,..., n)+0 (1) , (73) 

where Sj_j+i = 2pi ■ j9j+i and Sn,n+i = SnA- As usual fi is the scale introduced to maintain 
the dimension of the integral in 4 — 2e dimensions. 

When the implementation passes both non-trivial tests we have checked all coefficients 
proportional to infrared and coUinear divergent master integrals. This leaves only the 3- 
leg off-shell triangle coefficients and the 4-leg off-shell box coefficients unchecked. Note 
that for the 4-gluon and 5-gluon ordered amplitudes there are no finite master integral 
contributions and hence all coefficients are checked by looking at the divergent parts. The 
6-gluon amplitude gets a contribution from the finite 3-leg off-shell triangle master integral. 
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FIG. 6: The relative error for 100,000 ordered 4-gluon amplitude for the (+ H ) helicity choice. 

The horizontal axis is the log- 10 of the relative error of Eq. ()75p . the vertical axis is the number 
of events in arbitrary linear units. The left figure is the contribution, the middle figure is the 
contribution and the right figure is the finite part. 

To compare with the analytic results we generate 100,000 flat phase space events for the 
2 {n — 2) gluon scattering using RAMBO 33|. The center-of-mass collision energy is 
V^. The events are required to have the following cuts on the outgoing gluons: a cut on 
the transverse energy, Et > 0.01 x \/S, a maximum rapidity, rj < 3 and a separation cut, 
AR > 0.4. 

The evaluation time for 10,000 events is: for a 2 — 2 gluon ordered helicity amplitude 9 
seconds, for a 2 — > 3 gluon ordered helicity amplitude 35 seconds and for a for a 2 — > 4 gluon 
ordered helicity amplitude 107 seconds. Note that using the integration-by-parts method of 
ref. js] the evaluation time for 10,000 events would be approximately 90,000 second. This 
means the unitarity method of section 2 improves the evaluation of the six-gluon amplitudes 
by a factor of approximately 900, almost 3 orders of magnitude. The six-gluon evaluation is 
only three times slower than the five gluon evaluation and eleven times slower than the four 
gluon amplitude. This can be understood by counting the number of coefficients needed to 
evaluate the scattering amplitude. The number of coefficients multiplying a non-zero master 
integral for the n-gluon scattering ordered amplitude is 




— n 



24 ^ 



2n^ + lln - 34) 



(74) 



where the first term is the number of box coefficients, the second term the number of trian- 
gle coefficients and the third term gives the number of self energy coefficients. Finally, the 
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FIG. 7: Same as fig. [6l but the ordered 5-gluon amplitude for the (+ + 



-) hehcity choice. 



last term subtracts the n external gluon bubble master integrals (because the corresponding 
master integral is zero). That is, the number of 4-gluon coefficients is six, the number of 
5-gluon coefficients is twenty and the number of 6-gluon coefficients forty-four. The compu- 
tational time roughly follows the number of coefficients to be calculated for the scattering. 
This is a very different scaling law than the n\ growth of a straightforward Feynman diagram 
expansion (as is the case for the integration-by-part algorithm of ref. sl). The "factoriza- 
tion" of the diagrams into the tree-level blobs explains this large difference and causes the 
large gain in speed using the unitarity method compared to more conventional methods, 
replacing a n\ factorial growth by a ri^ power growth. The advantage of unitarity method 
with respect to a a Feynman-diagram based approach becomes greater as the number of 
external particles grows (provided the tree- level amplitudes are known). 

The numerical comparisons for the cut constructible part of m4(+H ), m5(+H ) 

and m6(+ H ) ordered helicity amplitudes are summarized in figs. (I6|7f8|) The 

degree of agreement is quantified by the expression 



S = log 



10 



m 



(1) 



unitarity "^analytic 



(1) 



m 



(1) 



(75) 



analytic 

In the 3 figures the 100,000 events are compared and binned in the quantity 5* for each of the 
3 contributions: e~^, and finite. As can be seen the majority of the events agree with a 
relative precision of 10"^ or better. This is more than sufficient for TEVATRON, LHC and 



We also compared all other helicity combinations with the known results in the literature, leading to 
similar results. 
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FIG. 8: Same as fig. [6l but the ordered 6-gluon amplitude for the (+ + 



-) hehcity choice. 



ILC applications. However, for the and finite parts a small portion of the events have 
a worse agreement. This is related to the fact that in these cases the amplitude receives a 
contribution from the bubble coefficient. In particular it is related to the calculation of the 
triangle subtraction term needed in the calculation of the bubble coefficient. For a small 
fraction of the phase space points there is a numerical instability in the matrix inversion 
needed to calculate the coefficients of the triangle spurious term in Eq. fl50l) . To understand 
this better we simply rotate the two basis vectors of the trivial space of the triangle in the 
following manner 



= n J* + i rig 



nz 



n, — in: 



2 5 



(76) 



such that 

n+ ■ 72+ = 72- ■ = 0; 72+ ■ r2„ = 2 
Within this basis the spurious term becomes 



C = c'q + c[s+ + + C^s\ + C4S_ + CgS?. + CgS^ , 



(77) 



(78) 



where s+ = I ■ n+ and s_ = I ■ ri-. By choosing 7 loop momenta {li}l=i "we get the set of 
equations 

EM^4 = E"^^A?, (79) 



mi 



i=l 



where rrii = c{li) and Aij = (l, s+(/j), s^(/i), si^(/j), s_(/j), s?.(/i), s"L(/i)). We note that the 
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matrix Aij is a double vanderMonde-matrix 43j, i.e. a matrix of the form 



1 ai al al A (3^ 
1 aa ai ai /^a /5| /3| 



(80) 



^ 1 07 /3| /9f J 

which is known to give an unstable inverse. No numerical procedure is known to stabilize the 
calculation of the inverse. In the figures we simply used the inverse anyway to see if these 
numerical issues would appear. We see that for a small fraction of the events they indeed 
appear in the case of the and finite contributions. On the other hand it is important 
to realize we have an infinite set of equations (i.e. an infinite set of loop momenta fulfilling 
the triple cut unitarity condition) with only 7 coefficients to determine. This means we can 
achieve arbitrary precision using a x^-type fitting procedure, at the cost of more computer 
time. We did not pursue this method for this paper. 

The final issue is the presence of Gram determinants in the box, triangle and bubble 
coefficients. We need to be aware of the Gram determinants when performing phase space 
integrals. We identify from the procedure in the last section two separate mechanisms of 
generating Gram determinant type of denominator factors in the coefficient. 

The first mechanism is straightforward. The solutions of the unitarity constraints for the 
box by quadruple cut gives a loop momentum which is proportional to the inverse of the 
square root of the 3-particle Gram determinant as can be seen in Eqs. ( !25ll27j) 

1 



(81) 



Because the maximum rank of the box tensor integral is 4, we can get at most terms of 
order A(/ci, /c2, k^)~'^. More precisely 

B C D E 



d. 



ijkl 



A + 



+ 



+ 



+ 



(82) 



A(A:i, k2, k^y/^ A{ki, k2, ks) A{ki, k2, k.fl^ A{h, k2, k,y ' 
with ki = {qj — Qi), k2 = {qu — Qj) and k-s = {qi — qu)- Similar, for the triangle and bubble 
coefficients we get respectively up to A{ki^ k2)~^^^ and A(ki)~^ terms. 

The second source of Gram determinants is more subtle and is generated for 5-particle 
or higher scattering amplitudes. It can happen that e.g. d^{li234) causing an instability 
at a specific phase space point. It should be treated with care when performing a numerical 
phase space integration. 
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IV. SUMMARY AND OUTLOOK 



We have presented a formulation of the 4-dimensional unitarity cut method in a physical 



language using the van Neerven-Vermaseren basis of ref. [30|]. This basis is also used to reduce 
5- (or higher) point tensor integrals to 4-point tensor integrals. When applied to 4- (or lower) 
point tensor integrals the decomposition of the loop momentum generates extra basis vectors 
spanning a "trivial" space. This trivial space is orthogonal to the external momenta of the 
loop integral. When applying unitarity cuts, the loop momentum dependence in the trivial 



space generate so-called spurious terms 



15| . These spurious terms integrate to zero when 



considering an individual cut diagram. However, when combining the double, triple and 
quadruple cuts care has to be taken not to double count. These spurious terms play an 
important role in the subtraction schemes needed to avoid the double counting problem. 
The procedure outlined in this paper allows us to solve the unitarity constraints without 
resorting to the explicit 4-dimensional spinor formalisms used in analytic calculations. This 
makes the method equally applicable to processes with (complex) masses for the internal 
lines. 

When applying the 4-dimensional unitarity cuts to the amplitude, both the master in- 
tegral coefficients and the spurious terms are calculable in terms of factorized products 
of tree-level amplitudes. The cut lines cause two of the external momenta of the tree- 
level amplitude to be complex. Existing leading-order generators, such as VECBOS 44 1 



and NJETS 



45 



46| . can be upgraded to allow for two complex external momenta with- 



out much effort. Once these upgraded tree- level generators are interfaced to the numerical 
program based on method described in this paper, they can be converted to one-loop gen- 
erators without any additional calculations. This will allow us to automate the calculation 
of the cut-constructible part of the one-loop amplitudes for many important background 
processes at the TEVATRON and LHC such as PP ^ 4 or 5 jets, PP ^ W +3, or 4 jets , 
PP ^t + t + l,2or3 jets , PP — >t + t + 6 + & + 0orl jets (massless bottom quark) etc. . 

As an example, we took the analytically known 4-, 5- and 6-gluon ordered tree-level 
helicity amplitudes and used them to evaluate the cut-constructible part of the 4-, 5- and 
6-gluon ordered one-loop helicity amplitudes. Especially the 6-gluon amplitudes are very 
complex, demonstrating the power of the method in the paper. The computational time 
required to evaluate the amplitudes is fast enough for serious TEVATRON, LHC and ILC 
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applications. The numerical instabilities in certain phase space points due to additional 
linear dependences between the particle momenta (the "Gram determinant instabilities") 
are easily identified in the master integral coefficients. The degree of the instabilities are the 
same as one gets from analytic unitarity calculations. 

The ultimate goal is to construct a NLO parton level generator. To reach that point, two 
outstanding challenges remain. The first one is the completion of the scattering amplitude, 
i.e. an automated calculation of the rational part. The second, far more difficult, challenge 
are the phase space integrations of both the one-loop amplitudes and the bremsstrahlung 
contributions. 
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